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Abstract 

We present a comprehensive analysis of hadron production in e~^e~ colhsions at 
different center-of-mass energies in the framework of the statistical model of the 
hadron resonance gas. The model is formulated in the canonical ensemble with ex- 
act conservation of all relevant quantum numbers. The parameters of the underlying 
model were determined using a fit to the average multiplicities of the latest mea- 
surements at ^/s = 10, 29-35, 91 and 130-200 GeV. The results demonstrate that, 
within the accuracy of the experiments, none of the data sets is satisfactorily de- 
scribed with this approach, calling into question the notion that particle production 
in e~^e~ collisions is thermal in origin. 



1 Introduction 



The analysis of hadron yields measured in central heavy ion collisions from AGS up to 
RHIC energies has shown [T1[2|3|^f5|6|7f8|9l[T0l[TT] that hadron multiplicities can be de- 
scribed very well with a hadro-chemical equilibrium approach which is governed by the 
chemical freeze-out temperature T, baryo-chemical potential /if,, and the fireball volume 
V; for a recent review see [I2]. The natural question arising here is whether this statis- 
tical behavior is a unique feature of high energy nucleus-nucleus collisions or whether it 
is also applicable in elementary collisions like, e.g., e~^e~ . Previous publications indicated 
that indeed hadron production in e~^e~ collisions at 14-43 GeV p^l4lll5] and 91 GeV 
[T31IT6] can be well described within a thermal model provided that local quantum num- 
ber conservation is properly implemented. The main result of these investigations was 
that the temperature values deduced are almost constant near T = 160 MeV and that 
the volume increases with energy, while strangeness is undersaturated. These results were 
taken, together with the results for nucleus-nucleus collisions where a similar temperature 
is reached at high energies, as evidence for the interpretation that the thermodynamical 



state is not reached by dynamical equilibration among constituents but rather is a generic 
fingerprint of hadronization [T71IT8] or a feature of the excited QCD vacuum [19]. Alterna- 
tively, it was argued in [2U] that the quark-hadron phase transition drives the equilibration 
dynamically for nucleus-nucleus collisions. Equilibration in e~^e~ collisions is not easy to 
explain in this latter approach. 

In our new analysis we are using the latest multiplicity measurements summarized and 
published by the Particle Data Group (PDG) [21]. Since the aim is a precision calculation 
for a small system we employ a fully canonical form of the statistical model [T2|22] con- 
serving baryon number N, charge Q, strangeness S, charmness C, and bottomness B. In 
the present analysis charmed and bottom hadrons are relevant only via their feed-down 
contributions to the yields of the lighter hadron species (for more details see below). To 
reach within the model a precision comparable to that of the data from the LEP collider (a 
few percent), we have performed computations including quantum statistics (see below). 



2 The model 



The canonical statistical model we will base our investigations on is described in 
[T^13fl4lll6f22] . Here we present a short summary with emphasis on the way the quantum 
number conservation is implemented. Most hadronic events in high energy e"'"e~ annihila- 
tions are two-jet events, originating from quark-antiquark pairs of the five lightest fiavors. 
Since we would like to address the issue of overall equilibration in these systems, we are 
using a scheme in which each jet is treated as a fireball with vanishing quantum numbers 
as fixed by the entrance channel; we call this the "no fiavor" scheme. It is clear at this 
point that hadrons from jets with heavy quarks (c and b) will be greatly underestimated 
by the model because of the large Boltzmann suppression factors. In this approach the 
issue of equilibration is effectively addressed only for hadrons with light quarks {u, d, s). 

It is important to recognize that the measured yields of these hadrons contain the con- 
tribution from the e"'"e~ annihilation events into cc and bb. Heavy-quark production is 
indeed significant and is very precisely measured, in particular at the Z° mass (-y/i=91.2 
GeV), where the measurements are very well described by the standard model [21]. Hence, 
heavy-quark production is manifestly non-thermal in origin. We therefore consider three 
scenarios: i) we fit the data as measured, ii) we subtract from the yields of hadrons car- 
rying light quarks the contribution originating from charm and bottom decays based on 
available data for charmed and bottom hadron production (and their branching ratios) at 
91 GeV. iii) we perform a fit to the data at 91 GeV in a 5-fiavor (u, d, s, c, b) approach. 
In this case e"'"e~ — > qq events are treated in as 2-jet initial state where each jet carries 
the relevant quantum numbers. The fractions of the quark fiavors in hadronic events [21] 
are thus external input values, unrelated with the thermal model (see also Table II in 
ref. [23]). We will treat the heavy quark sector in detail in a forthcoming publication [2^ . 

Scenarios i) and ii) differ from that used e.g. in ref. [23], where a 5-flavor scheme is 
used throughout. A significant difference between our approach and that of [23] is in the 
treatment of the volume entering the statistical model calculations, as explained below. 
In our case the volume is the hadronization volume of each jet, meaning that each jet 
hadronizes separately, as in deep inelastic scattering. The yields calculated for each jet 



are then added to compare with data. In ref. [23] the volume is that of both jets together. 
Because of different canonical corrections (see below) this leads to significant differences 
in the final results. While we consider our jet hadronization picture as more suited to 
describe two-jet events within a statistical frameworlj^. we acknowledge that the approach 
of hadronization in one volume [23] cannot be excluded at present. 

The appropriate tool to deal in a statistical mechanics framework with a system where 
all quantum numbers are zero is the canonical partition function with exact conservation 
of the N, Q, S, C and B quantum numbers 
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y^^) j=i k , . 
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where Xj is a five component vector = {Nj,Qj, Sj,Cj, Bj) containing the quantum 
numbers of hadron j, and = {(pjy, (f)Q, (ps, (pC: (ps) contains the parameters of the sym- 
metry group [U (1)]^. The sum over k is for the phase space cells for each hadron, f] = 1/T, 



^j,k = yPfe + rrij. In this expression, each (px corresponds to the conservation of the corre- 
sponding quantum number X; Nb and Np count the total number of boson and fermion 
states, respectively. 

The usual procedure is to solve Eq. [T]in the Boltzmann approximation, i.e. ln(l±x)^^ ~ x 
To perform the calculation with full quantum statistics we employ the series expansiorT^ 

ln(l-x)-i = f:^ (2) 
k=i 

in Eq. [1] leading to 

^(^) = y d'cP e^^^ exp{5: ^je""^^ + ^ zte-^'^^^} (3) 

with 
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where rrij is the particle mass and gj is the spin-isospin degeneracy factor. The index j 
runs over all particles species in the hadronic gas and b runs only over bosons. Here and 
in the following we use units with h = c = 1. The complete derivation of the partition 
function (for quantum statistics) will be the subject of a separate publication [2S]- The 
differences between calculations with Boltzmann and quantum statistics are presented in 
Table [H 

The set of hadron species used here was updated compared to [TT] with the most recent 
information available from the PDG ED and includes all resonances listed therajj. We 



^ Three and more jet events are neglected here as they are in other approaches too. 

^ We perform calculations with quantum statistics for bosons only. Due to the higher masses, 

for fermions the Boltzmann approximation is sufficiently accurate. 

^ We have included the a meson, as in ref. [27] 



have also, as in ref. [TT], included the resonance width explicitelj0- 



The integral representation of the partition function in Eq. ([3]) is not convenient for nu- 
merical analysis as the integrand is a strongly oscillating function. Thus we have applied 
an expansion into series of Bessel functions [28ll29f30] to obtain a result that is free of 
oscillations. We obtain the multiplicity (jij) per jet for particle species j by introduc- 
ing a fugacity parameter \j which multiplies the particle partition function and by 
differentiating 



d\nZ 



(5) 



One has to take into account that (rij) is the yield resulting from one jet. The multiplicity 
for the whole event is then the sum over the two jets. 

Table 1 

Comparison of particle yields obtained with our code and with the THERMUS code [31] for 
both initial production (prior to strong decays) and for the final values (after strong and elec- 
tromagnetic decays). We show yields (sum of particles and antiparticle yields and for the 2 jets) 
calculated with our code with quantum statistics (QS) and with Boltzmann statistics (BS) for 
the parameter set T=158 MeV, V=30 fm^ and 7^=0.80. 



particle 


this work 


THERMUS 


final 


initial 


initial 


final 


QS 


BS 


QS 


BS 


TT+ 


18.54 


18.07 


5.12 


4.74 


4.72 


14.28 




11.08 


10.73 


2.97 


2.68 


2.68 


8.29 


K+ 


2.028 


2.016 


0.954 


0.945 


0.940 


1.89 


KO 


1.952 


1.939 


0.938 


0.928 


0.924 


1.84 


ri 


1.090 


1.087 


0.472 


0.468 


0.472 


0.890 


/(770) 


1.12 


1.12 


0.756 


0.753 


0.756 


1.044 


K*°(892) 


0.597 


0.595 


0.447 


0.445 


0.442 


0.570 


P 


0.974 


0.966 


0.232 


0.230 


0.229 


0.668 


(j)il020) 


0.131 


0.131 


0.128 


0.128 


0.128 


0.132 


A 


0.364 


0.361 


0.0706 


0.0703 


0.0677 


0.239 


S+(1385) 


0.0393 


0.0389 


0.0313 


0.0310 


0.0300 


0.0316 




0.0202 


0.0200 


0.0115 


0.0114 


0.0108 


0.0190 


H0(1530) 


0.00768 


0.00760 


0.00738 


0.00731 


0.00679 


0.00679 


n 


0.00125 


0.00123 


0.00125 


0.00123 


0.00117 


0.00117 



Before engaging in the present data analysis we have performed a comparison of results 
from the model described above with those obtained using the THERMUS code [21], which 
is available publicly. The results for a particular set of parameters are presented in Table [H 



^ We have not used the method applied in [13|14|15|16j to include the uncertainties in hadron 
mass, width and branching ratios as additional systematic errors. 



Here and in the following the quoted yields include the corresponding antiparticle states 
(i.e. the yield labelled vr"*" is actually the sum of the yields of vr"*" and vr") and are the sum 
over the two jets, as summarized and published by the PDG [2T]. The agreement between 
the particle yields obtained with both codes is generally very good, lending strong support 
also to the numerical implementation of our methods. Note that, in our calculations, we 
employ quantum statistics, whereas in the THERMUS code the Boltzmann approximation 
is used for the canonical ensemble. Neglecting quantum statistics causes an error of about 
3% for the final pion yield, i.e. larger than the uncertainty in the data. As expected, the 
effect on all other hadron yields is smaller. If we employ the Boltzmann approximation 
in our code we get, prior to strong decays, values for mesons and non-strange baryons 
which are in agreement with THERMUS results at the percent level, while for strange 
baryons our yields are systematically higher by about 2% compared to those obtained with 
THERMUS. Inspecting the yields after strong decays, one notices a discrepancy between 
results from our code and from THERMUS for most of the hadrons. Our larger yields are 
due to a more complete set of hadron species used in the calculations. Concerning the 
decays, we have ensured that, in our code, the decay tables are symmetrical for particles 
and antiparticles and the decay widths always add up to the total width, even if particular 
channels are not measured. As the branching ratios are not well known for some of the 
high-mass states, we have used the known BR of the nearest state with the same quantum 
numbers. 



3 The fit procedure 

The multiplicity calculation proceeds in two steps. First, a primary hadron yield, Nfl^, is 
calculated using and ([5]). A crucial assumption of the model is that the final yields of all 
particles are fixed at a common temperature, the chemical decoupling point. As a second 
step all resonances in the gas which are unstable against strong and electromagnetic decays 
are allowed to decay into lighter stable hadrons, using appropriate branching ratios (B) 
and multiplicities (M) for the decay j ^ h published by the PDG [21]. The abundances 
in the final state are thus determined by 

Nh = iVf + E ■ Bij h)M{j ^ h) (6) 
j 

where the sum runs over all hadron species. 

In the resonance gas model the results are determined by the basic thermal parameters, 
temperature T and volume V (the volume corresponding to one jet). Following the ap- 
proach of ref . [13] , we also introduce an additional parameter 7^ into the partition function 
to account for a possible deviation of strange particle yields from their chemical equilib- 
rium values. If a hadron contains Ug strange valence quarks, its production is reduced by 
a factor 7"". This parameter is also applied to neutral mesons such as 77, 77', (j), u, /2(1270) 
and /2(1525) according to the fraction of ss content in the meson itself. The relevant 
fraction is determined using mixing formulae quoted in |21j . 

For the fit procedure we use the complete set of all measured yields of hadrons carrying 
light quarks. The second scenario, namely subtracting the contribution originating from 
charm and bottom decays, is investigated only at y^=91 GeV as a case study, since the 



measurements needed to allow the subtraction are complete only at this energy. The 
fit is performed by minimizing 

^^^^iNT^ (7, 

h 

as a function of the three parameters T, V and 7^, taking account of the experimental 
uncertainties ah. 



4 Results 



The resulting best fit to the data at the energy y/s=91 GeV is shown in Fig. [H Shown are 
the three cases of the fit, the non-flavor approach without and with the subtraction of the 
contribution from heavy quarks (see below for the magnitude of this contribution) and the 
5-flavor approach. We first note the overall behavior of the data, namely an approximately 
exponential decrease of particle yield with increasing particle mass. Such a behavior is 
expected in the hadron resonance gas model due to the Boltzmann factors, thus indicating 
the presence of statistical features of hadron production in elementary collisions. 

The quantitative description of the data with the statistical model is, however, rather 
poor and certainly no improvement is visible for the case of subtracting charm and bottom 
contributions. The poor fit quality which is already visible in Fig. [T] becomes striking when, 
as in Fig. [21 we show for the four energies the difference A (in units of the experimental 
error) between the experimental data and the statistical model calculations for the best 
fit values. 

A summary of the fit parameters obtained for the various data sets is presented in Table |2l 
The errors of the thermal parameters listed in Table [2] are the statistical errors as extracted 
from a minimization done with MINUIT [32], interfaced with our code. For the discussion 
of the systematic errors see below. 

Table 2 

Values of fit parameters and values per degree of freedom for different energies. At 91 GeV 
the fit results are given for the three scenarios discussed in the text. 



V^[GeV] 


T[MeV] 


V[fm3] 


7s 


xVdof 


10 


166±1.7 


10±1.5 


0.80±0.02 


314/20 


29-35 


162±1.7 


16±1.4 


0.96±0.03 


108/17 


91 i) 


160±0.5 


26±1 


0.82±0.007 


595/28 


91 ii) 


168±0.5 


18±1 


0.66±0.01 


994/28 


91 iii) 


170±0.5 


16±1 


0.66±0.01 


499/28 


130-200 


154±2.8 


40±4.3 


0.82±0.03 


12/2 



Typical values per degree of freedom lie between 5 and 20 and discrepancies between 
single data points and fit values larger than 5 standard deviations are not rare. Further- 
more, there is no clear pattern observed: the fits are comparably poor for baryons and 
mesons, as well as for non-strange and strange hadrons. In particular, for all energies the 
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Figure 1. Comparison between the best fit thermal model calculations and experimental hadron 
multiplicities (sum of particles and antiparticle yields and for the 2 jets) for e~^e~ collisions 
at y^=91 GeV. The upper panel shows the fit in the non-flavor scheme of data including 
the feed-down contribution from heavy quarks, the middle panel is after subtraction of this 
contribution, the lower panel is for a 5-flavor scheme fit where the flavor abundancies are extra 
input parameters from data (see text). The best fit parameters are listed for each case. 

yields of </> mesons and of hyperons are poorly reproduced. Large deviations are seen also 
for kaons. 

The case of the fit after subtraction of the contribution of the decays of charmed and 
bottom hadrons (scenario ii), explored at ^/s=91 GeV, is characterized by a larger 
value compared to the overall fit (scenario i). The extracted fit parameters differ somewhat 
in the two cases. In particular, T is higher and V is smaller for scenario ii) compared to 
scenario i), following the (T,V) correlation discussed below. Scenario iii) gives the lowest 
value, although the value is still far from that of a good fit. 
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Figure 2. Difference (in units of experimental error) between experimental data and thermal 
model fits at four energies. For ^5=91 GeV with open squares we show the results of the fit to 
data after subtraction of the heavy quark contribution (scenario ii), while the triangles are for 
the fit for the 5-flavor scheme (scenario iii). 

It should be mentioned that, while the general agreement among the four LEP experiments 
is excellent, the measured yields of S* hyperons differ by more than 70%. Excluding the S* 
hyperons would cause an increase of 7^ at ^/s=91 GeV, because the S* yields calculated 
in the model are overestimated (see Fig. [1]). This would slightly improve the situation for 
the S and A multiplicities which are higher than those predicted by the model, but with 
only a marginal improvement of the values, as is discussed below. 

In Table [3] we show for selected hadron species the calculations in the no flavor scheme 
and the experimental data after the subtraction of the charm and beauty contribution. 
The relative magnitude of this contribution is also listed in Table [31 

A difficulty for the determination of fit parameters is visible if one inspects the contour 
lines as shown in Fig. [3]in (T,V) plane for fits at 91 GeV. One notices in this figure a strong 
anticorrelation between the fit parameters which is also present in the (T,75) space (not 
shown). Closer inspection reveals, in addition, a series of local minima which indicates 
the difficulty in the determination of the fit parameters. Such local minima are typical 
for poor fits and imply that the true uncertainties in the fit parameters are likely much 
larger than the values obtained from the standard fit procedure [32] employed here. 



Despite these caveats about fit quality and uncertainties it is noteworthy that the tem- 



Table 3 

Calculated yields (corresponding to the fit of the subtracted data at 91 GeV, T=168 MeV, V=18 
fm^ and 7s=0.66) for selected hadron species for the non-flavor scheme and for the 5-flavor one 
assuming vanishing or fractional baryon and charge quantum numbers. In the third column we 
show the percentage of the yields arising from the contribution of charm and bottom events, as 
calculated from the experimental data (complemented by calculations using scenario iii) for He, 
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S+(1385) 


0.0337 


0.0232 


3.0 
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0.0231 


10.6 
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0.00056 


0.00140 


12.8 



perature parameters obtained from most data sets are close to 160 MeV and nearly in- 
dependent of energy, similar to results of previous investigations. In contrast, the volume 
increases with the center of mass energy. The values obtained for the strangeness un- 
dersaturation parameter 7^ range between 0.96 and 0.70 and exhibit no clear trend with 
energy. For the 5-flavor case, at 91 GeV we obtain a volume of 18 fm^. This is smaller 
than the results reported in [23] by about a factor of 2 because, in our case, the 2 jets 
hadronize separately. 

In the following we will take a closer look at the 91 GeV fit, as the data set at that 
energy contains the largest number of measured hadron yields. We consider scenario i), 
that is, the fit to data without subtracting the charm and bottom decay contributions 
to calculations in the non-fiavor scheme. To check whether the high x^/dof values are 
caused by discrepancies for a few particular particles, we excluded 3 hadron species (the 
E*'s and 0) and repeated the fit. Naturally, the fit is better, but we still found a high 
xVdof=208/25 (for T=158 MeV, V=28 fm^ 7,=0.90). 

The data from LEP comprise very many particle species and their yields typically are 
measured with accuracies of a few percent. Clearly the precision of the experimental data 
set provides a stringent test of the statistical model. To provide a quantitative estimate at 
which accuracy level the statistical model breaks down and to allow comparison with the 
situation encountered for RHIC data we performed a fit at y/s = 91 GeV to those hadron 
yields which were used for the fit of central nucleus-nucleus collision data at y/sNN=200 
GeV [m, namely 7r+, 7r°, K+, K*+, K*°, p, A, S", fi", and (j). The resulting thermal 
parameters are: T=162±2 MeV V=24±2 fm^ and 7,=0.82±0.02, with a still poor xVdof 




Figure 3. contour lines in temperature and volume space for the overall fit (left panel) and 
after subtraction of charm and bottom decays (right panel) for ^/s=91 GeV for the no flavor 
scheme. The contours correspond to +20%, +50%, +100%. The best fit values are 

indicated by the crosses. 

of 278/8. Excluding the 0-meson from the fit yields T=160±2 MeV, V=24±2 fm^ and 
7s=0. 96+0.02, with x^/dof = 62/7. A more direct comparison between thermal fits for 
heavy ion and e^e~ data can be obtained by arbitrarily assigning the uncertainties of the 
RHIC data to be corresponding LEP yields (namely uncertainties of the order of 10%). 
Constraining, as for the RHIC data, the fit parameters to T and V, i.e. setting 7^ = 1 
yields then T = 158+2 MeV and V = 24+2 fm^, with xV^of = 49/9. A reasonable fit 
can only be obtained by also letting 7., vary freely, with resulting parameters T= 168+2 
MeV, V=18+2fm3 and 7,=0.80±0.02, with xV^of = 18/8. These exercises demonstrate 
that a statistical model description of e~^e~ data fails badly without the introduction of 
the non-equilibrium parameter 7^. Even using 7^ the statistical description breaks down 
completely at an accuracy level for the data better than 10%. 

Another noteworthy difference between fireballs in e~^e~ and nucleus-nucleus collisions 
is their energy content. For its determination we have computed the energy density e 
in the hadronic gas to yield the thermal energy content E^ = eV of the jet at chemical 
decoupling. For the e~^e~ case and the parameters reported in Table [2l this procedure 
leads to Ej = 3.23, 5.45, 8.6, 10.2 GeV at = 10, 29-35, 91, 130-200 GeV, respectively. 
Hence, the thermal energy within each jet is only a small fraction of \/s/2, (e.g. about 
19% at 91 GeV). Apparently, in the thermal interpretation of e^e" collisions, most of 
the cm. energy is not available for particle production. This is in strong contrast with 
results for nucleus-nucleus collisions. We have analyzed for this purpose central collision 
events for 20 and 40 GeV/nucleon Pb-Pb collisions [11], where it makes sense to consider 
data integrated over the full phase space. In these cases we find that the energy content 
in the fireball amounts to 61 % and 63 % of the total cm. energy at 20 and 40 AGeV, 
respectively, implying that most of the total cm. energy in a nucleus-nucleus collision is 
thermal, with the remaining non-thermal fraction likely to be due to collective flow. We 
note that these differences are not consistent with the flnding that particle production 
in e+e~, pp and nucleus-nucleus collisions is universally governed by the available cm. 
energy [331IM1I35] . 



5 Summary and conclusions 



We analyzed the comprehensive set of measured yields of hadrons with light quarks {u, 
d, s) in e'^e~ collisions in a range from 10 GeV up to 200 GeV within an equilibrium 
thermodynamics picture. The calculations were performed in the framework of canonical 
partition functions to conserve explicitely baryon number, strangeness, charge, charm- 
ness and bottomness and made use of the hadron resonance gas description. Our results 
corroborate previous findings [T3lll4|15|[T6] that statistical features are present in hadron 
production in e"'"e~ collisions. At ^/s=91 GeV, three scenarios were considered, namely 
fitting data without and with subtraction of the decay products of charmed and bottom 
hadrons in our no flavor approach (vanishing jet quantum numbers) and employing a 
5-flavor scheme where the q and q jets (with the relative abundancies of the qq events 
as extra parameters taken from measurements at LEP) are carrying the flavor quantum 
numbers. The two thermodynamical parameters temperature T and volume V, and a 
strangeness undersaturation factor 7^ were obtained from a minimization procedure. 
While we find, as in previous investigations [T3lll4lll5lfT6] . that the resulting temperature 
value is close to 160 MeV, independent of energ\l ^ [the overall description of the high- 
precision LEP data is rather poor, independent of whether heavy quark contributions are 
subtracted or not. The x^/dof values larger than 5 for all fits call into question the validity 
of the thermodynamical approach for these data. This conclusion still holds even if the 
analysis is restricted to the same set of hadrons which were analyzed in the context of 
thermal model fits to Au-Au collision data from the RHIC accelerator. 

The apparent statistical fingerprint visible in the LEP data and first observed in [13] breaks 
down at an accuracy level of about 10%. Even at that level the e~^e~ data cannot at all 
be described without the explicit assumption of strangeness undersaturation, implying 
that hadron production in e~^e~ originates from a state which is quite far from true 
thermodynamic equilibrium. This conclusion is further supported by the observation that 
the corresponding fireball volume contains only a small fraction of the overall cm. energy, 
implying that most of the cm. energy is not available for particle production. This is in 
strong contrast to the situation in nucleus-nucleus collisions. 

A striking feature observed for nuclear fireballs is the complete absence of strangeness 
suppression. In this context it would be very interesting to probe the equilibrium features 
observed in central nuclear collisions with increasing accuracy. Will non-equilibrium fea- 
tures there also be first revealed in the strangeness sector? To clarify these issues will be 
a challenge for future measurements. 
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